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An efficient density matrix renormalization group (DMRG) algorithm is presented and applied 
to Y-junctions, systems with three arms of n sites that meet at a central site. The accuracy is 
comparable to DMRG of chains. As in chains, new sites are always bonded to the most recently 
added sites and the superblock Hamiltonian contains only new or once renormalized operators. 

Junctions of up to A = 3n + 1 500 sites are studied with antiferromagnetic (AF) Heisenberg 

exchange J between nearest-neighbor spins S or electron transfer t between nearest neighbors in half- 
filled Hubbard models. Exchange or electron transfer is exclusively between sites in two sublattices 
with Na ^ Nb- The ground state (GS) and spin densities pr =< S^ > at site r are quite different 
for junctions with S — 1/2, 1, 3/2 and 2. The GS has finite total spin Sg — 2S{S) for even (odd) 

N and for Mg = Sg in the Sg spin manifold, pr > 0(< 0) at sites of the larger (smaller) sublattice. 

S = 1/2 junctions have delocalized states and decreasing spin densities with increasing N. S = 1 
junctions have four localized Sz = 1/2 states at the end of each arm and centered on the junction, 
consistent with localized states in S' = 1 chains with finite Haldane gap. The GS of S' = 3/2 or 2 
junctions of up to 500 spins is a spin density wave (SDW) with increased amplitude at the ends of 
arms or near the junction. Quantum fluctuations completely suppress AF order in S' = 1/2 or 1 
junctions, as well as in half-filled Hubbard junctions, but reduce rather than suppress AF order in 
aS = 3/2 or 2 junctions. 

PACS numbers: 75.10. Pq, 05.10.-a, 75.10.-b 


I. INTRODUCTION 

The transport and magnetic properties of a system 
with a junction of three wires have been a frontier area 
of research. Y junctions such as 3-terminal Josephson 
devices [T] or carbon nanotubes [2] have been studied 
experimentally. Understanding quantum effect in three 
terminal junctions is important for potential applications 
as rectifiers m [3], switches and logic gate devices [4]. 
Recently these systems have also been studied theoreti¬ 
cally. [SHU] Interesting predictions include a low energy 
chiral fixed point with an asymmetric current flow in 
a spinless fermionic system m and negative density 
reflection at the junction of Bose liquid of ultra-cold 
atoms la. Theoretical studies have been mainly based 
on field theoretical approaches mils]. 

Exact numerical results are limited to very small junc¬ 
tions and tend to be inconclusive, especially with respect 
to quantum many body effects. Numerical techniques 
such as density matrix renormalization group (DMRG) 
give excellent results for one-dimensional (ID) systems 
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m- At first sight, however, DMRG appears to be far 
less accurate for structures with three terminals in which 
a long bond is repeatedly renormalized. Guo and White 
(GW) introduced [17] a new DMRG algorithm, summa¬ 
rized in Section II, for Y junctions with spin S at every 
site. We present in this paper a modified DMRG algo¬ 
rithm, quite distinct from that of GW, for Y-j unctions 
with N = 3n-\-l sites and three equal ID arms of n sites. 
The accuracy and efficiency of the modified algorithm is 
comparable to DMRG in ID chains, and we have studied 
Y junctions of up to 500 spins. We note that tensor-tree 
networks m are a general approach to many-body 
systems with a tree structure such as Y junctions, 
dendrimers or Bethe lattices. Stilbenoid dendrimers are 
a recent quantum chemical application [19] based on 
molecular units with many degrees of freedom. Tree net¬ 
works based on different units call for diverse algorithms. 

We consider Y junctions with antiferromagnetic 
Heisenberg exchange J between sites with spin S. The 
Hamiltonian for the junction in Fig. [^is 

Hg - J ^ ^ Srp • . (1) 

{rr') 

The sum is restricted to adjacent sites and J = 1 is a con¬ 
venient unit of energy. We discuss systems with S = 1/2, 
1, 3/2 or 2 and also study fermionic junctions that cor¬ 
respond to half-filled Hubbard models with N electrons 
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FIG. 1. Y junction of iV = 10 sites with three equal arms 
of n = 3 sites. The numbering for the left (unprimed), up 
(primed) and down blocks (double primed) is used in the 
DMRG algorithm; the numbering in red is used for spin den¬ 
sities. 


and N sites, 

Hp = —t {cLp^ap'o- + ^p'a^pcr) ^ ^pa^pjd^P^^pa • 

{pp')cr P 

( 2 ) 

Electron transfer t (set to 1 to define the energy scale) 
is limited to adjacent sites the number operator is 
Up, and U > 0 is repulsion for two electrons at a site. 
The Hiickel or tight-binding limit of U = 0 is readily 
solved exactly and provides a direct check of accuracy. 
Trimethylenemethane, C 4 ^Hq^ is a Y junction with 
four TT-electrons, four C atoms and has a stable triplet 
GS [20]. Bipartite Hiickel models with many singly 
occupied 7r-orbitals are a design principle for high spin 
hydrocarbons [21]. The atomic limit U » t reduces 
Hp dit half filling to Hs with S = 1/2 and J = /U. 

Both Hs and Hp conserve total spin Sp and its com¬ 
ponent Sz- By convention, we choose the Zeeman com¬ 
ponent Sz = Sp when the ground state (GS) has spin 
Sg > 0. The spin density at site r is the GS expectation 
value 


Pr = (S^r) ( 3 ) 

The sum over sites r returns 5';^ > 0, but individual pr 
may be positive or negative. Y junctions are bipartite: 
All exchange J or electron transfer t is between sites 
that form two sublattices, A and B, here with Na Np 
sites. The GS of Hs has Sq = S\Na — Nb\^ which 
alternates between Sg = S and 2S for odd and even 
respectively. Sites in the larger sublattice have positive 
pr, those in the smaller sublattice have negative pr. 


The paper is organized as follows. The modified al¬ 
gorithm for Y junctions with equal arms is presented 
and tested in Section |TT| including both infinite and fi¬ 
nite DMRG algorithms. Its accuracy is fully compara¬ 
ble to DMRG for chains. As in chains, new sites are 
always bonded to the most recently added sites and the 
superblock Hamiltonian contains only new or once renor¬ 
malized operato rs. The algorithm is applied to Y junc¬ 
tions in Section HI, first to fermionic and S' = 1/2 junc¬ 
tions, then to S = 1 junctions and finally to S = 3/2 
and 2 junctions. We focus on spin densities and size 
dependence. Localized states in S = 1 junction are in 
excellent agreement with the valence bond solid (VBS) 
model of Affleck, Kennedy, Lieb and Tasaki (AKLT) [22] . 
There is a localized Sz = 1/2 state at the end of each 
arm and one centered on the junction. The localization 
length ^ = 6.25 in arms is close to the chain result of 
White and Huse [23] while the length O < ^ indicates 
greater localization at the junction. The GS of S = 3/2 
or 2 junctions up to 500 spins are unexpectedly differ¬ 
ent, however: They are spin density waves (SDWs) with 
increased amplitude at the ends of arms and near the 
junction. Antiferromagnetic (AF) order is possible in 
systems with 5 'g' > 0 and found in Y junctions of S' > 1 
spins. Quantum fluctuations suppress AF completely in 
S = 1/2 or 1 junctions with long arms, but only partially 
in S > 1 junctions. We briefly mention in the Discussion 
generalizations of the algorithm to other junctions. 


II. MODIFIED DMRG ALGORITHM 

The computational effort for one eigenstate in con¬ 
ventional DMRG for ID chains with open boundary 
conditions goes as 0{Nm^) where N is the number of 
sites and m is the number of states per block for a given 
accuracy [H] [25]. The reason why is as follows: The 
number of arithmetic operations to obtain all eigenvalues 
of an L X 1/ matrix is 0{L^); so the number of operations 
for one eigenvalue goes as L?. In DMRG the matrix size 
is L = 16m^ for fermions and L = (25' + l)^m^ for spin 
5. In either case, goes as and a system of size N 
requires DMRG steps of N/2. This estimate excludes 
the construction and diagonalization of density matrices 
which are 0{m^). The greatest cost is the superblock 
diagonalization that goes as O(m^) for one eigenstate. 

Conventional DMRG for Y junctions scales as O(m^) 
and, as shown in Fig. of ref. 17, involves a long bond 
whose operators are renormalized many times. GW m 
cite previous DMRG applications to Y junctions and 
present a more efficient scheme for junctions with three 
equal arms that meet at a point, as in Fig. or at 
the vertices of an equilateral triangle. Singular value 
decomposition is used to obtain the density matrix of a 
single arm, and is faster as it requires O(m^) operations 
instead of O(m^). The method has a large truncation 
of the density matrix when two arms are combined 
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into a single block. We avoid truncation below since 
arms are never combined. The slow 0{m^) step is the 
GS eigenvector of the superblock matrix, which goes 
as L = 0{{2S l)m^) for three arms. The modified 

algorithm has improved accuracy for smaller m and 
makes it possible to treat Y junctions of Y ^ 500 sites. 

DMRG algorithms of ID chains add two sites per 
step between the left and right blocks. [24l [25] The new 
sites are always bonded to the most recently added 
sites. The superblock Hamiltonian of chains only 
contains new operators or once renormalized operators, 
a desirable featured that we retain for Y junctions. In 
some algorithms the number of newly added site in the 
superblock can vary from one [26], two [IT] or four |5] 
depending on the accuracy requirements of the systems. 
Here the superblock grows by three sites. 

The modified algorithm is shown schematically in Fig. 

Sites enclosed in loops define systems that consist of 
an arm plus a new site. The environment contains all 
sites in the other two arms. The key point is that the 
system at one step becomes an arm in the next step, 
thereby avoiding having to combine two arms into one 
block. The system block has m x r degrees of freedom, 
m for the arm and r for the new site, with r = 4 for 
fermionic junctions and r = 25' + 1 for sites with spin 
5. The relevant dimension is not mr of the system 
block, however, because the density matrix is block 
diagonal in sectors with different Ms values. The time 
needed to diagonalize the density matrix is negligible 
in sectors of dimension mr/{NS + 1). We obtain all 
density matrix eigenvectors \i) of the system block by 
block diagonalizing it into different Ms sectors. 

As shown in Fig. we start with a superblock of four 
sites, Y = 3n + 1 = 4. The notation (2,1,1) refers to an 
arm plus the central site and two other arms, respectively. 
The second step corresponds to (3,2,2) and Y = 7, the 
third to (4, 3, 3) and Y = 10, and so on. The n + 1 sites 
of arm plus central site at step n become the arm at step 
n + 1. We find the GS eigenvector of the superblock, 
starting with Y = 4, and expand \^) in the basis of the 
system (arm plus site) and the environment (two arms), 

\'>P) = '^Cik\i)\k). (4) 

ik 

The basis vectors \k) are direct products of basis states 
of two arms of the system of the previous step. The 
total number of sites, Y = 3n + 1, increases by three at 
each step. The reduced density matrix of the system has 
elements 


diagonalization, we take the m eigenvectors of p with the 
largest eigenvalues as elements of an Mxm matrix p'. 
The effective Hamiltonian and operators in the truncated 
mxm basis are renormalized according to 


= {p')^Op 


where {p'Y is the transpose matrix, and O and H are 
the operators and Hamiltonian of system block. The 
superblock eigenvalue calculation is the slow step that 
scales as O(m^), although conservation of total 5^ 
reduces the dimension to less than (25 + l)m^. The GS 
then yields the reduced density matrix pij of the system 
for the junction in which each arm is one site longer. 
Since operators of the system block are renormalized 
only once, similar to ID chains, we expect similar 
truncation errors in Y junctions. 

The following steps and Fig. describe the infinite 
DMRG algorithm for Y junctions with equal arms: 

(a) Start with four sites, the superblock in Fig. 

(b) Find the GS eigenvalue and eigenvector. 

(c) Gonstruct the density matrix of the system block, 
shown in Figs[^, and[^ for 4,7 and 10 site su¬ 
perblocks, respectively. Diagonalize it to get the 
eigenvectors corresponding to the m largest eigen¬ 
values. 



Pij=Y.^*ikCjk. (5) 

k 

The sum is over the environmental degrees of freedom. 
We suppose p to be a matrix of dimension M. After 


FIG. 2. Schematic representation of the infinite DMRG al¬ 
gorithm for Y junctions with equal arms. At each step, the 
loop encloses the system and the superblock contains a new 
site shown as an open dot, and three arms. 
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(d) Renormalize the operators and Hamiltonian for the 
system blocks using Eq. 

(e) Construct the Hamiltonian of the superblock as 
shown Fig. and 

(f) Repeat the process from b to e until the desired 
size N = 3n 1 is reached. 

Finite DMRG is required to obtain accurate spin den¬ 
sities, correlation functions and other GS properties. The 
conventional finite algorithm for ID chains has two new 
sites and sweeps through two arms of the same chain 
[24l|25]. The algorithm for Y junctions has one new site 
and sweeps through two arms while keeping the third arm 
constant. The procedure is shown schematically in Fig. 
and summarized below. Three to four DMRG sweeps are 
typically sufficient for converged energies. Finite DMRG 
is particularly important for junctions with S = 3/2 and 
2 sites. 

(a) Start with the superblock with equal arms as shown 
Fig. ^ taken from the infinite algorithm calcu¬ 
lation. Select two arms A and B for sweeping 
through. 

(b) Find the GS eigenvector of the superblock. The 
new system block ‘A’ is the old block A plus a new 
site, shown as the open dot in Fig. Wp- Block ‘A’ 
has one site added to arm A and removed from 
arm B at every step while arm G remains same. 
Gonstruct the reduced density matrix of the system 
block ‘A’. 

(c) Renormalize the Hamiltonian and operators of the 
system block ‘A’. 



FIG. 3. Schematic representation of finite DMRG steps. Sites 
of the system block are enclosed in the loop; remaining sites 
are in the environmental block.The new site is the open dot. 
A, B and C arm refer to the three different arms. 


TABLE I. The m dependence of the GS energy per site Sq of 
Y junctions with 64 sites, equal arms, and U = 0 in Eq. |2lor 
S = 1/2, J = 1 in Eq. [g 


m 

o 

II 

S = 1/2 

20 

-1.23336877828 

-0.43915791387 

40 

-1.25809972580 

-0.43915891861 

60 

-1.25826370430 

-0.43915892503 

80 

-1.25838125000 

-0.43915892523 

100 

-1.25842968750 

-0.43915892525 

Exact 

-1.25848468281 

- 


(d) Gonstruct the superblock Hamiltonian with new 
blocks ‘A’, ‘B’ and ‘G’ shown in Fig. [^. 

(e) Repeat the steps b to d until block ‘B’ has only one 
site as shown in Fig. 

(f) Now take ‘B’ as the system block, block ‘A’ and ‘G’ 
as the environment. As shown in Fig. |^, add a 
site to block ‘B’ and remove one from block ‘A’. 

(g) Repeat steps b to d for system block ‘B’ until block 
‘A’ has only one site as shown in Fig. |^. 

(h) In the next step, ‘A’ becomes the system block. 
One site is added in ‘A’ and removed from block ‘B’. 
Steps b to d are repeated until Fig. is reached. 

(i) Now take blocks ‘A’ and ‘G’ while keeping block ‘B’ 
constant and repeat the cycle a to h that starts and 
ends with equal arms in Fig. ii. Finally, take block 
‘B’ and ‘G’ while keeping the block ‘A’ constant. 
Repeat the cycle from a to h that starts and ends 
with Fig. 

(j) One cycle of finite DMRG is the whole process from 
a to i. 

Next we discuss the accuracy and efficiency of the 
algorithm. The U = 0 limit of Eq. is a Hiickel or 
tight-binding model of non-interacting electrons on N 
sites that can readily be solved exactly. As an example, 
we took a half-filled band of Y = 64 sites and calculated 
the GS energy per site eo as a function of m, the 
dimension of the system block in the truncated basis. 
Table [l] shows good convergence by m ^ 60 for this 
fermionic system of about 4^ degrees of freedom, or 
some 4^^ per arm. DMRG of non-interacting electrons 
often converges the most slowly due to GS degeneracy or 
to higher entanglement entropy [26l[28]. The Y junction 
of S' = 1/2 spins in Eq. is the U ^ t limit with 2^ 
spin degrees of freedom whose GS energy per site is not 
known exactly. As shown in Table [Ij m ^ 20 is sufficient 
for So of junctions with 21 spins per arm. 

As additional tests of the algorithm, we consider the 
total energy E{m) of 64-site Y junctions with J = 1 in 
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FIG. 4. Total GS energy F^(IOO) — F^(m) as a function of m for 
64-site Y junctions with equal arms, J = 1, and the indicated 
S per site in Eq. 


Eq. and S = 1/2, 1, 3/2 and 2 as a function of m. 
The truncation errors P{m) = 1 — Xj on keeping m 
eigenvalues of the density matrix are listed in Table [ll| 
for S = 3/2 and 2. 

Since the exact GS is not known, we follow the evolu¬ 
tion of AE{m) = E{mo) — E{m) where ttiq = 100 is the 
nominally the converged value. Excellent convergence is 
achieved in Fig. [^by m ^ 70, with AE of the order of 
10-10 for S = 1/2,10-^ for S = 1, and 10-^ for S = 3/2 
or 2. Increasing m to 130 lowers AF^/F^(100) by 5 x 10“^ 
for S' = 2. The P{m) change in Table is also small. 
By contrast, the GW algorithm [17] for AE with S = 1 
reaches only 10-^ around m = 140 in Fig. 1^ of ref. 17. 
The present algorithm is well suited for Y junctions, both 
because as in ID chains operators are renormalized only 
once and because the procedure in Fig. [^increases the 
number of sites smoothly without ever having to combine 
two arms. 


III. LOCALIZED STATES AND 
ANTIFERROMAGNETIC ORDER 

We apply the modified DMRG algorithm to Y junc¬ 
tions with equal arms, either half-filled junctions in Eq. 


TABLE II. Truncation errors P(m) of 64-site Y junction of 
S = 3/2 and 2 as a function of m. 


m 

P(m),S = Z/2 

P(m),S = 2 

64 

1.2 X 10“® 

7.2 X 10““ 

80 

7.6 X 10“^° 

3.0 X 10“® 

100 

1.7 X 10“^° 

1.5 X 10“® 

130 

5.2 X 10““ 

5.6 X 10“'^ 


or Heisenberg junctions with spin S at every site in 
Eq.[T] Unless otherwise stated, the results are based on 
m = 100 and 5 — 10 sweeps of finite DMRG. We discuss 
junctions of A = 3n + 1 sites, distinguish between 
odd and even A, and study the size dependence. The 
algorithm is applicable to junctions of A ^ 500 sites. 
We focus on AF order in Heisenberg junctions with 
S > 1 and on localized states of junctions with integer 
S'. 

As mentioned in the Introduction, Y junctions are bi¬ 
partite, with different number of sites Na Nb in sub¬ 
lattices A and B. We take Na > Nb and have 

A + 1 

Na = —-— = Nb + 1 (odd A, even n) 

N 

Na ~ ~2 ^ ^ ~ ^ ^ (even A, odd n). 

The junction is in sublattice A for odd A, in sublattice 
B for even A. The Neel state \AE) has spins PS at all 
sites in sublattices A and B, respectively, and is the SDW 
with the largest possible amplitude; \AE) is exact in the 
limit of classical spins, S' ^ oc. Quantum fluctuations in 
Eq. [^strongly reduce AF order for S > 1 and suppress 
it altogether for S = 1/2 or 1. Nevertheless, \AE) gives 
the correct spin, Sq = 2S for even A and S for odd A, 
and also accounts for the sign of the GS spin densities in 
Eq. with pr > ^ for r in sublattice A and Pr ^ 0 for ^ 
in sublattice B. 


A. Fermionic and S = 1/2 junctions 

The Hiickel junction has U = 0 in Eq. [^and Na — Nb 
nonbonding orbitals with energy 5 = 0 and nodes at 
all sites in sublattice B. The nonbonding orbitals are 
easily found analytically. The half-filled junction has A 
electrons, A sites and spin a in nonbonding orbitals. 
The GS for odd A has Sq = 1/2 and pr = 1/{2Na) 
at sites in sublattice A, = 0 at sites in sublattice B. 
The triplet GS for even A has So = Pr = 1/A^ at 
sites in Na and = 0 at sites in Nb- Since 5'^ = 1 for 
arbitrarily large (even) A, the Hiickel densities at sites 
in sublattice A decrease as 2/(A + 2). 

Increasing U > 0 in the half-filled junction does not 
change So but induces negative < 0 at Nb sites and 
increases p^ > 0 at Na sites. The sum over |p^| increases 
with U as localized spins are formed due to electron 
correlations. The spin densities are no longer equal, 
however, as seen in Fig. [^at U = 4t, the bandwidth of 
the ID Hiickel or tight-binding model. The Heisenberg 
model with S = 1/2 in Eq. has the largest positive 
and negative spin densities. 

The Heisenberg junction with even A has pj < 0 at the 
junction and spin densities that go as 1/A. The S = 1/2 
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FIG. 5. Spin densities pr along any two arms of 64-site Y 
junctions with U = 0 or 4t in Eq. or S' = 1/2 in Eq. 
The junction is at r = 22 with pj < 0. 


junction of = 202 spins also has a triplet GS and a spin 
density distribution similar to the Heisenberg model in 
Fig-IU Longer arms lead to smaller spin densities: pes = 
—0.0929 for the junction at r = 68; pe? = 0.0974 and 
p66 = —0.0817 at the first and second neighbors of the 
junction; pi = 0.0307 and p 2 = —0.0198 at the first two 
sites of arms. As expected, quantum fluctuations entirely 
suppress AF order in the infinite S = 1/2 junction. We 
note that spin densities increase along the arms at odd r 
and become more negative at even r. We will later And a 
different pattern in S' = 3/2 junctions in which quantum 
fluctuations are not as dominant. 


B. S = 1 junctions 

Haldane [29] predicted finite energy gaps A(S) in infi¬ 
nite Heisenberg spin chains with integer S and nearest 
neighbor J > 0. Experimental realizations of S = 1 
chains have confirmed a gap that DMRG evaluates [23| 
as A(l)/J = 0.4105. The valence bond solid (VBS) 
picture of AKLT [22] has been widely applied to S = 1 
chains, and we do likewise for S = 1 junctions. S = 1 
chains with open boundary conditions have a localized 
state with = 1/2 and localization length ^ = 6.03 [23] 
at each end. GW m obtained four localized Sz = 1/2 
states in a Y junction with N = 181, one at the end of 
each arm and one centered on the junction. 

Figure [^ shows the spin densities in one arm of a Y 
junction of A" = 202 sites with S = 1. As expected, the 
GS has Sg = 2, the junction at r = 68 has pj < 0, and 
the total Sz of either localized state is 1/2. The spin 
densities in the first 15 sites of an arm and 14 sites from 
the junction are listed in Table [ni| The spin densities of 
localized states are conventionally taken as proportional 



FIG. 6. Spin densities in one arm of a Y junction of A = 202 
sites, S' = 1, as a function of r with r = 1 at the first site and 
at r = 68 at the junction. 


to [saisn] 

Pr cx {-lY~'^exp{-r/£,) (8) 

where r = 1 refers to the ends of chains. This ap¬ 
proximation neglects the difference between p^ for even 
and odd r that is clearly seen in Fig. and Table [IH| 
Any pair r, r + 2 defines a local localization length ^ 
^ = 2/{ln\pr \ — ln\pr^2\)‘ As seen in Fig. the > 0 
and pr <11 series have similar localization whose average 
is ^ = 6.25 for arms, = 5.81 for the junction, and the 
first few sites deviate from a simple exponential. White 
and Huse [23] obtained ^ = 6.03 for S = 1 chains with 
open boundary conditions; they did not consider positive 
and negative spin densities separately. GW m report 
similar localization at the junction and arms without 
going into detail, while we And slightly but distinctly 
smaller = 5.81. 

We consider next Y junctions of S' = 1 spins and odd 
A = 199. The GS is a triplet. So = Sz = 1, and the 
junction has pj > 0. Quite remarkably, the spin densities 
of the localized state are identical to a part per 10^ to 
the A = 202 values in Table IHII aside from a reversed 
sign around the junction. The localization lengths = 
6.25 and = 5.81 obtained for A = 202 are equally 
applicable to A = 199 within our numerical accuracy. 
Spin densities near the junction add to three localized 
states at the ends of arms for A = 202 and Sz = 2^ while 
they subtracts for A = 199 and Sz = 1- Identical \pr\ 
for A = 199 and 202 directly confirm that each localized 
state has Sz = 1/2. The GS of a Y junction of S' = 1 spins 
and long arms is 2^ = 16-fold degenerate and comprises 
a quintet, three triplets and two singlets. The quintet 
has A symmetry under Cs, the singlets transform as E, 
and the triplets as A and E. 
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TABLE III. Spin densities of a Y junction of N = 202 sites 
with S' = 1. Listed are the first 15 sites of an arm, the junction 
and up to 14 sites from the junction. 


spindensitypr 

Arm, r — 1 

Junction 

1 

0.5321 

-0.3044 

2 

-0.3209 

0.3530 

3 

0.3733 

-0.2515 

4 

-0.2652 

0.2459 

5 

0.2624 

-0.1886 

6 

-0.2000 

0.1737 

7 

0.1855 

-0.1383 

8 

-0.1469 

0.1234 

9 

0.1317 

-0.1004 

10 

-0.1068 

0.0880 

11 

0.0939 

-0.0726 

12 

-0.0773 

0.0629 

13 

0.0671 

-0.0629 

14 

-0.0558 

0.0450 

15 

0.0480 

-0.0377 


C. Junctions with S > 1 

The Haldane gap of the infinite spin-2 chain is smaller, 
[3l] A(2) = 0.0886 ± 0.0018, about A(l)/5 and is less 
accurately known than A(l). The ends of S = 2 chains 
are expected to have localized S = 1 states with corre¬ 
spondingly larger Schollwock et al. [30] have discussed 
the S = 2 chain in detail using DMRG, quantum Monte 
Carlo and exact diagonalization methods; they interpret 
results in a VBS framework and report [301 Eg limited 


r-68 



FIG. 7. Spin densities from Fig. [^plotted as ln\pr\ vs. r 
at an arm and at the junction. Circles (squares) represent 
negative (positive) spin density. 


agreement. Chains of A = 270 spins with increasing 
m (to 180) cover a 25-fold spin density change. Fig. 
of ref. 28 shows ^ = 2/(ln\pr\ — ln\pr^ 2 \) to vary as 
^ 40 ± 10 up to r ^ 30 and to be almost constant, 

^ ^ 50, in the range 40 < r < 125. 

Smaller A(2) and localization ^ ^ 50 in S' = 2 chains 
indicate that Y junctions with longer arms are needed 
to study localized S = 1 states. Instead of localized 
states, however, and in sharp contrast to S = 1/2 or 1 
junctions, we find substantial AF order in both S = 2 
and S = 3/2 junctions as shown in Fig. [^for A = 448 
(left panel) and 298 (right panel). The junction is at 
r = 150 or 100, respectively. The spin density in the 
interior of arms oscillates between ±c at odd and even 
r. The amplitude increases at the junction and at the 
end of arms, in contrast to the spin densities of the 
S = 1/2 junction in Fig. whose magnitude decreases 
from the junction. The similarity in the behavior of 
the S = 3/2 and S = 2 junctions is noteworthy since 
the infinite S = 3/2 chain is gapless unlike the S = 2 
chain. We consider the main features together before 
pointing out differences between S = 3/2 and 2 junctions. 

By definition, Heisenberg exchange is between local¬ 
ized spins S at every site. The sum over p^ is the z 
component of spin in the given state. The sum over \pr\ 
normalized to NS is the fraction of unpaired spins; the 
Neel state \AF) with ±5' returns (AS')“^^^ |p^| = 1. 
We interpret SDW amplitudes c in Fig. [^in the interior 
of arms as AF order c/S that increases with S. The 
fraction of unpaired spins in S = 3/2 junctions is 0.293 
for A = 245 and 0.302 for A = 448; the fraction for 
S = 2 spins is 0.326 for A = 445 and 0.323 for A = 448. 

By contrast, the fraction is less than 0.1 in S' = 1 
junctions for A = 202 or 199 and clearly vanishes in 
the infinite junction since unpaired spins are in localized 
states. The fraction of unpaired spins also goes to zero 
in S = 1/2 junctions with increasing A as discussed 
earlier. 

We generalize Eq. |^for p^ to reflect the different be¬ 
havior of spin densities near the junction and ends of 
arms. We study positive and negative pr separately but 
use the same length for simplicity, ^ for arms and 
for the junction. The spin densities in one arm run from 
r = 1 to n, with A = 3n + l and the junction at r = n + 1. 
For even A (odd n), we analyze the spin densities accord¬ 
ing to 

p2r-i = c + aexp{-{2r - 1)/^) + ajexp{-{n + 2 - 2r)/^j) 

P 2 r = -{c + bexp{-2r/C) + bjexp{-{n + 1 - 2r)/0))- 

( 9 ) 

The first sum over odd sites has r = 1,2,.., (n + l)/2; 
the second sum runs to (n — l)/2 since there is one fewer 
even site. Both positive and negative spin densities for 
odd A (even n) have r up to n/2 in Eq. 
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FIG. 8. Spin densities pr in one arm of Y junctions with (a) 
N — 448 and (b) N — 298 spins S' = 3/2 and 2. Lines are 
Eq. [^with parameters in Table [In| 


Exponential contributions are limited to either end 
when n > ^ or ^j. The parameters are obtained as 
shown in Fig. for = 448 (n = 149) and S = 3/2. 
The SDW amplitude is c = 0.336, is significantly 
smaller than ^ and ±pr lead to nearly equal ^ or ^j. 
The junction has pj = —0.510 for N = 448 and 0.480 
for N = 445; the arms have equal spin densities within 
a percent or two. The N = 445 junction has essentially 
the same parameters in Eq. Nevertheless, we have 
Sg = 3 for AT = 448, Sq = 3/2 for N = 445. The 
difference is largely due to the junction and its first 
few neighbors. Both Y junctions with even or odd N 
support a SDW with equal pi > 0 at the end of all three 
arms. Figure shows exponential contributions for 
S = 2 junctions with N = 448 and 5 'g = 4. Since the 
slow decrease of spin densities in Fig. [^and the resulting 
^ ^ 32 are in the expected range for A(2) ^ A(l)/5, the 
SDW amplitude of S' = 2 junctions may not be entirely 
due to end effects. We leave open the behavior of much 
longer junctions. 


Table [IV| lists the parameters of Eq. for junctions 
with even N that generate the lines in Fig. essentially 
the same parameters hold for A/'T 1. The size dependence 
of the S = 3/2 junction has apparently saturated or 
almost saturated ai N ^ 450, but has not saturated 
for the S = 2 junction. We always find O < faster 
decrease of the SDW amplitude near the junction. 


The VBS picture has localized S = 1 states at ends of 
S = 2 chains or the arms in S' = 2 junctions. We find 
decreasing pi = 1.202 and 1.200 for N = 298 and 448 
junctions while Schollwock et al. ED] report 1.13 for a 
chain A/^ = 270 with a fixed S = 1 defect at the other 
end. The next site has p 2 rsj — 0.95. The first 5-10 
spin densities deviate significantly from Eq. and we 


TABLE IV. Eq. [^parameters for Y junctions of N sites with 
spin S 


S,N 

c 

e 

O 

a, aj 

b, bj 

3/2, 148 

0.433 

5.2 

5.7 

0.38, 0.25 

0.24, 0.20 

3/2, 298 

0.350 

11.2 

9.1 

0.30, 0.21 

0.24, 0.19 

3/2, 448 

0.336 

12.2 

7.9 

0.33, 0.22 

0.26, 0.27 

2, 298 

0.68 

12 

10 

0.43, 0.39 

0.33, 0.34 

2, 448 

0.54 

32 

20 

0.58, 0.50 

0.52, 0.47 


do not know how to identify a localized state. The 
S = 3/2 junction for N = 298 and 448 has decreasing 
pi = 0.781 and 0.780 that, perhaps coincidentally, is 
again slightly larger than S/2 = 0.75. Since the SDW 
amplitude is S in Neel state \AF), quantum fluctuations 
reduce AF order by 50% at the ends of arms and by 
more than 50% elsewhere. SDWs occur naturally in sys¬ 
tems whose GS has 5 'g > 0 and 25 'g + I degeneracy in Sz- 

To conclude this Subsection, we comment on S' = 3/2 
and 2 chains with open boundary conditions that were 
motivated (i) by the unexpected result that S = 2 
junctions do not follow VBS and {ii) to confirm quan¬ 
titative agreement with Schollwock et al. [30]. The 
GS of quantum chains with an even number of spins 
A^ is a singlet, Sg = 0. It is not degenerate, thereby 
excluding a SDW, but may have quasi-long-range order 
in the infinite chain. Delocalized states are expected in 
the gapless S = 3/2 chain. The gapped S = 2 chain 
may have localized S = 1 states at either end that 
become decoupled in the infinite chain. Two localized 
states lead to exponentially small gaps between the 
singlet GS, a triplet and a quintet, just as S' = 1 chains 
have an exponentially small gap to the lowest triplet 
[33] . Accordingly, we studied the quintet. So = 2, 
with the lowest energy of S = 2 chains and for compar- 


r-150 



FIG. 9. Spin densities from Fig. [^ for S = 3/2, N = 448 
plotted as ln{\pr \ —0.336) vs. r at an arm and at the junction. 
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r-150 



FIG. 10. Spin densities from Fig. for S' = 2, N = 448 
plotted as ln{\pr\ — 0.54) vs. r at an arm and at the junction. 



FIG. 11. Spin densities pr in half chains of (a) A/" = 150 
and (5) N = 300 spins S = 3/2 and 2 with antiferromagnetic 
Heisenberg exchange J between neighbors 


ison the lowest-energy triplet, So = 1, of S = 3/2 chains. 


Spin densities for open S = 3/2 and 2 chains of 


N = IbO and 300 spins are shown in Fig. 
middle, where they are zero by 
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up to the 
symmetry. In either 
case, the first few pr depend weakly on size, as found 
previously [30], and are almost the same as in = 450 
junctions with 150-site arms. End effects are similar 
in chains and junctions, and exponential fits over a 
limited range are possible aside from the first few spin 
densities. Symmetry about the middle of chains leads 
to linear pr around r = N/2 as shown in Fig. The 
fraction of unpaired spins is large: 0.183 and 0.121 for 
S' = 2, AT = 150 and 300; 0.090 and 0.067 for S = 3/2, 
N = IbO and 300. We infer that the spin densities are 
primarily due to end effects in these chains or junctions, 
in sharp contrast to the localized states in S = 1 chains 
or junctions. 


IV. DISCUSSION 

We have presented a modified DMRG algorithm for 
Y junctions in Section [11] and results in Section |III| for 
junctions up to 500 sites, mainly junctions in Eq. [^ 
with Heisenberg exchange J between spins S = 1/2, 1, 
3/2 or 2. Much longer chains of, say, 1000 sites greatly 
increase the computational effort at the finite DMRG 
step. The accuracy may not be lower, however, since 
the entanglement entropy [28] of the GS for dividing 
the junction into system and environmental blocks 
will increase only slightly. As already noted, we are 
considering large but finite junctions rather than the 
thermodynamic limit. That limit is better studied in 
chains since neither the junction nor the ends of arms 
should matter in junctions with infinitely long arms. 


The accuracy of the modified algorithm is fully equal 
to the DMRG accuracy for ID chains. Two arms are 
never combined into one and new sites are always 
bonded to the most recently added sites. As in chains, 
the superblock Hamiltonian contains only new or once 
renormalized operators. Infinite DMRG is accurate for 
S = 1/2 or 1 junctions, where finite DMRG makes 
minimal improvements, but finite DMRG significantly 
improves the results for S' = 3/2 or 2 junctions. Three 
or four sweeps of finite DMRG is sufficient for good 
energy convergence. We performed 5 — 10 sweeps for 
spin densities in order to confirm the different GS of 
S = 1 and 2 junctions. 

The modified algorithm for Y junctions of equal arms 
can be generalized to other systems, to be discussed 
elsewhere [35]. All generalizations are based on the 
schematic procedures in Eigs. andfor equal arms, {i) 
No change is required for more than three equal arms, 
although computational requirement increase as dis¬ 
cussed on Section [n] on going from chains to three arms. 
{ii) The algorithm performs well in preliminary tests of 
Y junctions with arms of different lengths n ^ n' ^ n" 
[35] , The infinite algorithm with equal arms is run until 
the longest arm n is reached. Einite DMRG is then done 
using blocks of different size to construct the superblock. 
{in) GW considered Y junctions with arms that meet at 
an equilateral triangle instead of a point m- Eor such 
systems, the modified infinite algorithm can again be 
used to generate the desired junction. In the beginning 
of finite DMRG, the superblock is constructed using 
blocks of different size; [35] two blocks have the same 
size, the third block has one fewer site, and the new site 
is added to the third block. The modified algorithm can 
also be generalized to {iv) Y junctions with different S] 
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we use four new sites with different S and three arms at 
every step m- A new site is added at the end of each 
arm and another one is added at the junction of these 
three arms in Fig. 1^ Since the size of the density matrix 
is L = {2S + l)^mMor adding one spin S', the size scales 
as (2S + 1)^ for adding four spins S. The procedure is 
efficient for small S but rapidly becomes more expensive 
for large S. 

Kekule diagrams are special cases of VB diagrams in 
which singlet pairing is limited to adjacent atoms or to 
adjacent S = 1/2 sites. The complete VB basis also 
has singlet pairing of more distant sites. There are Na 
Kekule diagrams for the S = 1/2 junction with odd N: 
the unpaired electron is at a site of the larger sublattice, 
which uniquely fixes singlet pairing, shown as lines to 
neighbors, of all remaining sites. Although the total 
number of VB diagrams is exponentially larger, Kekule 
diagrams are often a simple approximation for fermionic 
or spin-1/2 systems that provides semi-quantitative 
information; they have been extensively used in organic 
chemistry for well over a century. Their scope now 
extends to VBS in which S > 1 sites are represented 
[22l|30l|34] as two or more S = 1/2 spins that are singlet 
paired with neighbors. VBS has 2S lines to neighbors 


and no unpaired spins aside from chain ends. Unpaired 
spins and resonance among Kekule diagrams account 
naturally for localized Sz = 1/2 state in S' = 1 junctions, 
albeit without any reference to the Haldane gap. 

VBS wave functions are exact GS of special Hamil¬ 
tonians that Schollwock et al. [30] discuss and write 
explicitly for S = 2 chains. Heisenberg chains with finite 
A(S) for integer S have bilinear exchange JSr • 
between neighbors while Hyss contains terms up to 
{Sr • for S = 2. Schollwock and Jolicoeur |32| find 

that the Haldane phase described by VBS is strongly 
reduced in S = 2 chains compared to S = 1 chains. 
Our DMRG results for finite Y junctions are quite 
consistent with VBS for S' = 1 junctions where quantum 
fluctuations suppress AF order. DMRG for S = 2 
junctions does not follow VBS, however. The GS is 
instead a SDW with reduced but finite AF order. 
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